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Basal-Plane Magnetic Anisotropics of High-K rf-Wave Superconductors 
in a Mixed State: A Quasiclassical Approach 
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We study the basal-plane anisotropies of reversible magnetization and torque in a mixed 
state of layered d-wave superconductors based on the quasiclassical version of the BCS-Gor'kov 
theory. Both the longitudinal magnetization (Ml) and torque (r) show fourfold oscillations as 
a function of the field angle \. The relationship between the node position and the oscillatory 
patterns shown by Ml and r is clarified. It is also shown that the sign of the r(x)-oscillation 
does not change between H c i and H c 2, while the sign of the Ml (x)-oscillation changes. The 
newly obtained result for r indicates that the torque experiment can allow us to detect the 
in-plane anisotropies of H C 2 even in a material with strong fluctuations such as cuprate or 
organic superconductors, where the H C 2 itself cannot be determined experimentally. 
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1. Introduction 

It has been widely known that the effects of non- 
local electrodynamics in a clean type-II superconduc- 
tor cause several macroscopic anisotropy phenomena. 1 
For the uniaxial anisotropy in a layered superconductor, 
or the basal-plane anisotropy in a multicomponent su- 
perconductor, the local London 2 or the local Ginzburg- 
Landau 3 ' 4 (GL) approach provides a simple and conve- 
nient description of the phenomenon in terms of the sec- 
ond rank mass tensor. However, as pointed out by Ho- 
henberg and Werthamer, 5 the anisotropy in a cubic sys- 
tem cannot be explained within the framework of these 
local theories, and in such a case, we must recourse to 
theories incorporating nonlocal effects. 6-10 Among these 
nonlocal theories, the quasiclassical Eilenberger formal- 
■ ism 10-12 is advantageous for a quantitative description 
£h °f the phenomena. It can also be employed to interpo- 
late the results obtained using the London theory and the 
GL theory, valid at H <C H c i and H « H c2 , respectively, 
since both theories can be derived from the Eilenberger 
formalism. 

Recently, in the field of unconventional superconduc- 
tivity, there has been renewed interest in such a mixed- 
state anisotropy phenomenon. Typical examples are the 
field- angle-dependent oscillations of the specific heat 13 ' 14 
and thermal conductivity. 15 These phenomena can be 
used to clarify the node position of the gap function, 
through the fact that the Doppler-shifted quasiparticles 
cause the variation of these quantities depending on the 
angle between the applied field direction and the nodal 
direction. For more details, we refer the reader to the 
seminal theoretical work by Vekhter et al. and a more 
quantitative calculation by Miranovic et al. 17 

Magnetization is another fundamental quantity in the 
mixed state of a type-II superconductor. Thanks to the 
newly developed "shaking" technique 18 and the arrival 
of high-quality low-pinning samples, we can now obtain 
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reversible magnetizations in a broad region of the field- 
temperature phase diagram. Indeed there exist precise 
reversible magnetization measurements 19,20 for borocar- 
bides, showing clear fourfold basal-plane anisotropies. 
In our previous work, 21 we studied such a basal-plane 
anisotropy of the longitudinal magnetization Ml. Wc 
demonstrated that the experimental data of Ml for boro- 
carbides can be explained by considering both gap and 
Fermi surface anisotropies. However, results on the trans- 
verse component Mp of the reversible magnetization, or 
torque r, have not been discussed. 

Experimentally, the basal-plane magnetic torque in a 
mixed state was measured by Ishida et al. 22 for an un- 
twinned single crystal, YBa 2 Cu 3 07_5 (YBCO). They ob- 
served a fourfold oscillation, in addition to a rather large 
twofold oscillation presumably due to the Cu-0 chains 
specific to YBCO. Referring to the result based on a 
nonlocal GL theory, 9 they concluded that the observed 
fourfold anisotropy is consistent with the simple d x i_ y i- 
wave model. However, their interpretation strongly re- 
lies on the conjecture that the sign of the r-oscillation 
does not change between H c i and H C 2 23 In this con- 
nection, it is worth emphasizing that the sign of the 
Mi-oscillation does change between H c \ and £f C 2. 21,24 
Namely, the anisotropy of Ml in low fields is opposite 
to what is expected from the i? C 2-anisotropy. Thus, it is 
of importance to clarify the r-oscillation behavior in the 
whole field range from H c \ to H C 2- 

In this work we study the basal-plane longitudinal 
magnetization Ml and torque t in a mixed state of lay- 
ered d-wave superconductors within the quasiclassical 
Eilenberger formalism. The present paper is an exten- 
sion of our previous work 21 on the longitudinal magne- 
tization to the transverse component, and clarifies that 
the sign of r-oscillation does not change while the sign of 
Mi-oscillation changes. Then, it is shown that a simple 
interpretation connecting the presence or absence of sign 
changes to the anisotropy dependences of H c \ and H C 2 is 
available. In addition, we provide the details of our anal- 
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ysis that were not discussed in our previous work. We use 
a simple d x 2_ y 2-wave model with the Fermi surface being 
isotropic within the basal plane. Since no detailed theo- 
retical study on the basal-plane magnetic torque has been 
reported, it is meaningful to study the simple d-wave 
model. Of course, in an actual material, another source 
of complication 21 may be the basal-plane anisotropy of 
the Fermi surface. For YBCO, the simple d-wave model 
is appropriate for the qualitative study of the nonlocal 
effects, since the d x 2_ y 2-w&ve model can explain the ob- 
served configuration of the vortex lattice in high fields, 

1. e., a square lattice with the nearest neighbors along the 
nodal directions. 25 

The plan of the present paper is as follows. In § 2, we 
develop a method of solving the Eilenberger equations 
for high-K superconductors. In addition, we introduce 
the approximate solution used in our previous work, 21 
as well as the solution given by Pesch. 26 By comparing 
these three results, we show that the approximate solu- 
tion presented in this paper gives a good result, at least 
with respect to the thermodynamics. This section may 
be skipped by those who are not so interested in the the- 
oretical details. In § 3, we introduce a model describing 
the basal-plane anisotropy of a layered superconductor. 
Based on the model and using the approximate solution 
presented in § 2, we discuss the basal-plane anisotropics 
of the longitudinal magnetization and torque, and show 
that the sign of Mi-oscillation changes, while the sign of 
r-oscillation does not. At the end of § 3, we give a simple 
interpretation on our microscopic results. The summary 
is given in in § 4. 

2. Quasiclassical formalism 

We use the quasiclassical Eilenberger formalism 11, 12 
derived from the BCS-Gor'kov theory through an expan- 
sion in terms of (pfCo) 1 *C 1 ; where p-p is the Fermi 
momentum and £o = v-p/2irT c with v-p and T c being 
the Fermi velocity and the transition temperature at a 
zero field. We start with the following Eilenberger equa- 
tions 11 for the quasiclassical Green's functions /, p and 
9 (&b = ft = c = 1): 

f(s n ,p,r) = (2e n + iv-n)- 1 (2g(e n ,p,r)w p A(r)j 
P(e n ,p,r) = f*(s n ,-p,r), 



where 



(...) 



J FS dS(---)/\v\ 
J FS dS/\v\ 



(5) 



denotes the average over the Fermi surface with dS being 
the area element of the Fermi surface. 

2.1 Full solution 

In this subsection, we propose a new method of solving 
the Eilenberger equations for the vortex lattice states in 
extreme type-II superconductors with large GL param- 
eters, k » 1. Since this method is beyond either the 
approximate solution in our previous work 21 or Pesch's 
solution, 26 the solution obtained by this method is called 
a "full solution" in this paper. 

We use the Cartesian coordinates (x\, x 2 , Xs) consti- 
tuting the right-handed system, and set the magnetic 
induction B parallel to x 3 . We divide the microscopic 
field b(r) into a spatially uniform part B (i.e., induc- 
tion), and a periodic part 8b(r) with zero mean within 
the vortex lattice cell: 



b(r) = B + 5b(r). 



(6) 



Since the latter part is known to be of the order of 
0(1/ k 2 ) <C 1 and negligible in high-K superconductors 
except for H w H c \, we neglect 5b from now on. To be 
consistent with this assumption, we neglect the demag- 
netization effect. The set of primitive translation vectors 
m, n of the vortex lattice in this paper is given as 

m = (2n/i>)xi, (7) 

n = (2n(/v)xi + vx 2 , (8) 

where lengths are measured in units of tb = (2|e|B) -1 / 2 
above and in the following. The real constants C and v 
specify the form of vortex lattice. For a triangular vortex 
lattice, we set £ = 1/2 and v = (37T 2 ) 1 / 4 , while for a 
square vortex lattice, we set ( = 1/2 and v = y/n. In 
any case, the difference in vortex configurations has a 
negligible influence on the thermodynamics at the level 
of the mean field approximation. 

We fix the gauge of the vector potential as 

,(1) 



A(r) = —Bx 2 X\. 



(9) 



g(s n ,p,r) = \Jl- f(e n ,p,r)p(e n ,p,r). 

The gap function is expressed as A p (r) = w p A(r) where 
w p is the pairing function with relative momentum p of 
the Cooper pair, and A(r) is the pair potential with cen- 
ter of mass coordinate r. Here, v = vpp with p = p/\p\, 
II = — iV + 2|e|A with A being the vector potential, 
and e n — 27rT(n + l/2) denotes the fermionic Matsubara 
frequency. In the following, only the region e„ > is con- 
sidered. We set w p — 1 in an s-wave superconductor and 
w p = V2(Px — py) in a c?-wave superconductor. The pair 
potential should be determined self-consistently through 
the gap equation 

(ln(f ) + 2ttT ]T e,; 1 ) A = 2ttT £ (w* p f), (4) 



(2^he boundary condition for the pair potential A(r) is 
gauge dependent, and with the present gauge, A(r) can 

(3be expressed as the superposition of each Landau level 
function TpN- 



Mr) 



N 



A ^2 d N ip N , 

N=0 

oo 



(10) 



H N (x 2 +vm) c _i {x2+v 



where A = 1.764T C , C m = (^) 1/4 e- i7rt > 2 , and H N 
is the N-ih Hcrmite polynomial. On the other hand, the 
normal component of the quasiclassical Green's function, 



E n >0 
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g(e n ,p, r), is a periodic quantity of r, and it can be ex- 
panded in the lattice sums by Fourier transformation: 

g{e n ,p,r) = ^.g Q (e„,p) cxp(iQ • r), (12) 
Q 

9q(£u,p) = \g(e n ,p,r)exp(-iQ -r) , (13) 

L J s 

where [• • • ] s denotes the spatial average. The reciprocal 
lattice vector Q is given by 

Q = mvx\ + (— mC, + n)(2ir/v)x2, (14) 

with m and n being integers. Substituting eqs. (10) and 
(12) into eq. (1), we have an expression of /, 



/ = 2w p ^2 g Q dN 0Q, N > 

Q,N 

Q , N = A (2e„ + i- i ;.n)- 1 (e iQ -> A r) 



(15) 
(16) 

,N- 



The remaining task is to derive the expression of 
Using the parameter representation (2e n + iv ■ II) _1 = 
/ °° dpe-^+^^P , wc obtain 



(2e r . 



r 

ru-n)- 1 = / 

Jo 



dp e 



-2e n p -iv 1 v 2 p 1 /r\ 



x g-i/roirii e -ipti 2 n 2 



(17) 

where the operator identity e A+B = e~^ A ' B ^e A e B was 
used, which holds if [[A,B],A] = [[A,B\,B\ = 0. Then, 
by applying the identity e a9x /(x) = f(x + a) for a 
nonsingular function /, and introducing the notation 
s = (p/r B ){viXi + v 2 x 2 ), we have 



N 



A r 

Jo 



dp e 



-2e„p 



xe 2 



SiS2+iS!X2 



tp N (r - 



JQ(r-s) 



(18) 



and we can finally solve the Eilcnbcrgcr equation (1). 

To numerically obtain the self-consistent solutions for 
/, f', and g, we use the following procedure. First we give 
initial values for {d^} and {pq}. Next we use eqs. (16) 
and (18) to obtain the new /, and eqs. (2), (3) to obtain 
P and g. From g(e n ,p,r), gQ{s ni p) is obtained using 
cq. (13). Each {djv} is determined by the following gap 
equation projected onto each Landau level: 

C £„>0 n > U £„>0 

(19) 

Finally, wc return to cq. (1) and repeat the same proce- 
dure until self-consistency is achieved. 

Once we correctly determine /, p, and g, the physical 
quantities can be calculated. In particular, the thermo- 
dynamic quantities can be obtained from the free energy 
F/V through the thermodynamic relation. The expres- 



sion ofF/V = F/V - B 2 /8n is given as 
1 , 



F 
V 



N(0) 



W\ 



£n>0 



2ttT (A*/ + A/t 



f(2e n + iv ■ U)p + p{2e n + iv ■ IT)/ 
2(l+ 5 ) 

(g-l)(f(w p Ay + p(w p A)) 



(20) 



nTN(0) 

e n >0 . 



(.9 + 1) 



(21) 



Here, IVI" 1 = ln(T/T c ) + 2nTJ2 en>0 £ n X is the strength 
of the attractive interaction. To proceed to the last line, 
we assumed that /, p, and A satisfy the Eilenberger 
equations (eqs. (1), (2)) and the self-consistent equation 
(eq. (4)). 

2.2 Approximate solution 

To solve eq. (1) in a more numerically tractable man- 
ner, we introduce the approximate solution that was used 
in our previous work. 21 Hereafter, the solution obtained 
by the method presented in this subsection is called the 
"approximate solution" . 

Let us introduce an approximation which is the gen- 
eralization of the one used by Pesch, 26 

f(e„,p,r) w 2g(e n ,p,r)w p ((2e„ + it;-n)- 1 A(r)). 

(22) 

This approximation relies on the fact that the normal 
component of the quasiclassical Green's function g has 
a weaker spatial dependence than the anomalous com- 
ponent /, as pointed out by Brandt et al. 27 We em- 
phasize here that in contrast to Pesch's solution, we do 
not replace g(s n ,p,r) in eq. (22) by its spatial average 
[g{ £ n,P, r)] s , therefore, we can reproduce the first non- 
Gaussian term of the nonlocal GL free energy given in 
rcf. 28. As we show later, Pesch's solution can be ob- 
tained from the approximate solution by completely ne- 
glecting the spatial dependence of the (/-function and by 
restricting A to the lowest (N — 0) Landau level. Re- 
cently, Pesch's solution has been widely used in slightly 
different contexts. 29-33 

To proceed with our approximation, it is convenient to 
define a new function <!> in a similar manner as described 
in ref. 34 (^-parameterization): 

*(e„,p,r) - 2 Wp ((2e„ + ru-II)- 1 A(r)). (23) 

This function offers simple representations of /, p , and 
9- 

I = .9$, (24) 
P = .9^, (25) 

9 = 7Twm> (26) 

where &(e n ,p, r) = <&*(£„, —p, r). Repeating essentially 
the same argument as described in the previous subsec- 
tion, we obtain 



$ = 2w p ^2 0Q=o,n d N . 



(27) 



N 
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After some calculations, /3q=o,jv in eq. (27) can be trans- 
formed into the form used in ref. 21, 



0Q=O,N 

a N (p) 



/>0O 

A / dpe- 2 ^a N (p), (28) 
Jo 

E H N (x 2 + vm - Re A) 



After some algebra, this quantity is transformed into 
P(e n ,p) = 4\w p d \ 2 Al 



/ dpi / 

Jo Jo 



dp 2 e 



-2 £ „( P1 +p 2 )-(pi+P2) 2 /4T^ 



m— — oo 

xe-d A ' 2 



VWw. 



-A 2 )/4 -(x 2 +vm-\) 2 /2-ivmx 



f-OC 

= 4\w p d \ 2 (A a T B ) 2 / dppe 
Jo 



p 2 /4-2e„r B p 



(34) 



(29) 

where ajy depends on p through A = (v 2 +ivi)p/rB- The 
numerical procedure for the above approximate method 
was explained in our previous work. 21 

Now we show that the approximate solutions of 
cqs. (24)-(26) can reproduce the nonlocal GL free en- 
ergy derived in ref. 28 up to the quartic term. We start 
from eqs. (24)-(26). Expanding /, f\ and g in terms of 
A, we have 

1 



ft = $t_l($t)2$ + .. 

.9 = 1- + . . . . 



(30) 



Then, we substitute these equations into eq. (21) and 
obtain 



F 
V 



N(0) 



A* ( —A 



ttT 



E n >0 



4 



+C.C. 



up to 0(|A| 4 ). This expression of F/V is equivalent to 
that derived in ref. 28 using diagrammatic perturbative 
calculation. This fact justifies the approximate solution 
near H c2 . 

2.3 Pesch's solution 

It is instructive to see how the approximate solution 
presented above is connected with Pesch's solution. 26 
Pesch's solution can be derived from the approximate 
solution by imposing two additional conditions: (i) com- 
pletely neglect the spatial dependence of the (^-function, 
and (ii) restrict A to the lowest (N = 0) Landau level. 
Below, we briefly explain how to derive Pesch's solu- 
tion from the approximate solution in a two-dimensional 
case. Following Pesch (see eq. (7) of ref. 26), we rewrite 
eq. (26) in the form <7 2 (1 H-<I><I>^) = 1. Keeping the above- 
mentioned two conditions and taking the spatial average 
of both sides of the equation, we obtain (cf. eq. (27) of 
ref. 30) 

= (l + P(e n ,p) 



\g(e„,p,r)] s 
P(£ n ,p) 



A\wpdo\ 2 Aj 



dp 2 e 



dpie 



-2e„p 2 



-2e„pi 



(32) 
(p)) 
«o(-P))]33) 



a { 



Here, p = (pi+ p 2 ) I 'tb with tb = rs/vp- Using Dawson's 
integral 35 W(i£) = erfc(£) with £ = 2e„rs, we finally 
have 

P(£ „, p)= (^^)^2( 1 _^ (i0 V (35) 

The above expression is equivalent to the result 30 re- 
formulated for a two-dimensional superconductor. The 
method of obtaining the self-consistent solution using 
Pesch's approximation is described in ref. 30. 

2.4 Validity of approximate solution 

In this subsection, we compare the approximate solu- 
tion in § 2.2 with the full solution in § 2.1 and Pesch's 
solution in § 2.3. Calculations are performed in two- 
dimensional cases with isotropic Fermi surfaces, where 
the area element of the Fermi surface is given by dS — 
p-pdep with (p being the azimuthal angle of v. First, we 
briefly explain our numerics. To perform the p-integral 
in the calculation of /3q : at, we use the simple trapezoidal 
rule. In our calculation, we use the cutoff for reciprocal 
-lattice vector Q in eq. (14) as — /Q max < m,n < /Q max ; 
(3]thus the total number of reciprocal lattice vectors equals 
(2/Q max +l) 2 . For the Fermi surface average, we use Simp- 
son's rule with an equally spaced mesh. For the spatial 
average, we use the simple trapezoidal rule, and we set 
iV max = 12. 

Now we show the result for the free energy density 
F/V = F/V — B 2 /8n. The magnetization can be de- 
rived from F/V through the thermodynamic relation 
M = —Vb{F/V). Figure 1 shows the field depen- 
dence of F/V for an s-wave superconductor and a d- 
wave superconductor. Here, the field B is measured in 
units of the two-dimensional orbital limiting field H 2 ^ = 
0.561(|e|/7r)/27r£ 2 . A triangular vortex lattice was as- 
sumed in the s-wave case, while a square vortex lattice 
was used in the d-wave case. We plot the results of the 
different methods discussed above: the approximate solu- 
tion; the full solutions with /Q max = 1 and 2; and Pesch's 
solution. 

As was discussed by Brandt et al. 27 the convergence 
by changing the cutoff Iq max for the reciprocal lattice 
vectors is quite excellent when we calculate the full so- 
lution of the Eilenberger equations. That is, the higher 
Fourier components of the g-function, <?q^o, do not give 
significant contributions to the thermodynamic quanti- 
ties. Indeed, we can see in Fig. 1 that data with /Q max = 1 
and 2 are almost the same. Moreover, data corresponding 
to the approximate solution and the full solution with 
Ig max = 2 cannot be distinguished on this scale. This 
means that the approximate solution provides a fairly 





Fig. 1. Field dependences of the free energy at T/T c = 0.5 (a) in 
an s-wave superconductor, and (b) in a d-wave superconductor. 
The squares, dash-dotted line, solid line, and dotted line repre- 
sent the approximate solution, the full solution with iQ max = 1, 
the full solution with iQ max = 2, and Pesch's solution, respec- 
tively. Here, = o.561(|e|/7r)/27rfg (?o = v f /2ttT c ) is the 
two-dimensional orbital limiting field. 

good result, at least with respect to the thermodynam- 
ics. On the other hand, Pesch's solution seems to overes- 
timate the free energy. This is because Pesch's solution 
integrates out the vortex degrees of freedom in advance, 
so that the system tends to approach the Meissner state. 

Next, we study the spatially resolved quantities. Figure 
2(a) shows the field dependence of {d^} at T/T c = 0.5 
in an s-wave superconductor. In the figure, a slight dif- 
ference appears between the full solution and the ap- 
proximate solution. For the lowest Landau level compo- 
nent do, we cannot see the essential difference between 
the full solution and the approximate solution. On the 
other hand, a non-negligible difference can be seen in 
the field dependence of higher Landau level components. 
We show, in Fig. 2(b), the corresponding data in a d- 
wave superconductor. Since we assume a square vortex 
lattice in this case, {d^v} with N = 0, 4, 8, . . . grows with 
decreasing field, 36,37 and this can be clearly seen in the 
figure. In both cases the approximate solution seems to 
slightly overestimate the amplitude of the higher Landau 
level components. Physically, the higher Landau level 
components describe the fine spatial structure of A, and 
these are relevant to the estimation of fine structure in 
the local density of states, which is closely related to 



Fig. 2. Field dependence of each {rfjv} at T/T c = 0.5 (a) in an 
s-wave superconductor, and (b) in a d-wave superconductor. The 
empty circles with the solid line represent the approximate solu- 
tion, and the filled circles with the dashed line represent the full 
solution with iQ max = 2. 

the differential conductance detected in the STM exper- 
iment. 38 These physics are, however, beyond the scope 
of this work. It is worth noting that the overall behavior 
of each {dw} calculated by the present scheme shows a 
good correspondence to the result obtained by Watanabe 
et al. 39 (see Fig. 2 therein), who solved the Eilenberger 
equations by a different method. A small difference in the 
^-component between our result and that in ref. 39 in 
the d-wave case is considered to stem from the difference 
in the configuration of the vortex lattice. We assumed a 
square lattice with the nearest neighbor oriented along 
the nodal direction, while in ref. 39, the nearest neighbor 
was oriented along the antinodal direction. 

The results in this section can be summarized as fol- 
lows. The approximate solution in § 2.2 gives relatively 
good results as long as we consider the thermodynamics 
in the mixed state. Therefore we can use this approxi- 
mation for the analysis of the basal-plane magnetization 
and torque, which we present in the next section. 

3. Reversible Magnetization and Torque 

3. 1 Model for layered superconductor 

In this subsection, we introduce the basal-plane 
anisotropy of a layered superconductor into our treat- 
ment. First of all, we discuss the anisotropic model that 
takes into account a layered structure of the materials, 
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Fig. 3. Coordinates (x,y,z) and crystal axes (a,b,c). Induction 
B is rotated from the fe-axis by angle \- 





(a) 



(b) 



Fig. 4. Vortex lattice structure seen in (a) real space coordinates, 
and (b) distorted coordinates. 



such as cuprate superconductors. An example of such a 
procedure 41,42 is the parameterization of the Fermi sur- 
face as a distorted cylinder. We consider a sample ar- 
rangement as shown in Fig. 3. We start with the uniaxial 
dispersion relation 



6p = 2m** +p2y ' ) ~ icos ( PzS )' 



(36) 



where s is the interlayer spacing and t is the inter- 
layer hopping energy. As usual, we define the uniaxial 
anisotropy parameter by 7 = vp/ts, where vf = pp/m. 
Assuming 7 S> 1 and pps 3> 1, the Fermi surface is pa- 
rameterized as 



p = p F p + p z z, 



(37) 



where p — cos(ip)x + sm(ip)y, —ir < ip < ir is the az- 
imuthal angle, and — n/s < p z < ir/s. Then, the Fermi 
velocity v = V p e p can be expressed as 



v = vf (p + 7 1 sin(p z s)z) . 



(38) 



and the area element dS divided by \v\ is given by 
dS/\v\ = md(pdp z . In this paper, we treat a simple 
d x 2_j,2-wave model, and we set the pairing function 

w p = cos (2(ip + x)), (39) 

where \ is the field angle, as shown in Fig. 3. 

Next we define the anisotropic Landau levels suitable 
for the description of the in-plane vortex state. To con- 
nect the formulation in the previous section with the 
present sample arrangement (see Fig. 3), we should set 
x\ = z, x 2 = x and X3 = y. In the following argument, we 
still use the notation [x\, x 2 , x 3 ). We define the raising 
and lowering operators a± as 

a± = -ifniiifjna), (40) 



2rj 

where the real constant i] characterizes the anisotropy. 
If we choose r\ > 1, then the vortex state defined 
by the operators is compressed along the xi-direction 
and stretched along the ^-direction (see Fig. 4(a)). Of 
course, the isotropic Landau levels are defined by set- 
ting 77 = 1 . Now it will be convenient to introduce a new 
coordinate, r' = (x^jX^), as 



Tr 



( V ^x uV -^x 2 ). 



Then the raising-lowering operators d'± and the resultant 
vortex state, expressed in the distorted r'-frame, appear 
isotropic (see Fig. 4(b)). Therefore, a natural way to de- 
fine the anisotropic Landau levels ip(r) in the r-frame is 

by 

Vw(r) = Vjv(T r) = i! N {ri 1/2 x 1 ,T 1 - 1/2 x 2 ), (42) 

where is the isotropic Landau levels defined in 
eq. (11). This definition is the extension of the idea 
used to describe the in-plane vortex state within the 
anisotropic lowest Landau level. 43-45 We note here that 
our treatment contains higher Landau level components. 
The important point in the above definition of ?pN(r) 
is that if we work with the r'-system in our calculation, 
the total effect of the uniaxial anisotropy can be absorbed 
into the anisotropy of the Fermi velocity represented in 
the distorted r'-frame: 



V- 1/2 v 2 - 



(43) 



Thus, all the formulations presented in the previous sec- 
tion can be valid, as long as we set x\ = z and x 2 = x 
in the following numerical calculation. Although the pa- 
rameter 77 can be chosen arbitrarily in principle, the nat- 
ural choice for r\ is to set r\ = 7. This is because if we 
construct a linearized anisotropic GL equation 46 from 
the model eq. (38), then the Landau levels eq. (42) with 
77 = 7 are the solutions of it under the field parallel to 
the a- or b- axes. 

3.2 Longitudinal magnetization 

In general, magnetizations can be derived from F/V = 
F/V — B 2 /8tt through the thermodynamic relation M — 
—Vb{F/V). Then the longitudinal magnetization -Ml(|| 
H) is obtained by Ml = —-^{F/V), since the direc- 
tion of the applied field H can be approximated by that 
of induction B in high-K superconductors. For this com- 
ponent, however, Klein and Pottinger 47 obtained a more 
convenient formula that is the extension of the virial the- 
orem derived by Doria et al.: 48 



-AitM l = 



47r 2 A^(0) f 
B 



E 

e„>0 



9(/ f KA) + /KA)*) 



(41) 



-2e„(l-$), 



(44) 
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Fig. 5. Logarithmic field dependence of longitudinal magnetiza- 
tion Ml in a d-wave superconductor. We set T/T c = 0.5 and 
7 = 3. Here, -ff or b = 1.037(|e|/7r)/27r£g is the three-dimensional 
orbital limiting field in the isotropic case. The inset shows the 



linear field dependence of M£ v 



M 



B II antinodc 



)/2- 



relatively wide range of field, as well as the GL behavior 
Ml oc H C 2 — B just below H c2 . This is because the GL 
theory is derived 49 from the quasiclassical Eilenberger 
formalism through an expansion about A, whereas the 
London theory is derived 6 using a phase-only (London) 
approximation. Thus we confirm that our analysis cap- 
tures the essential behavior in the mixed state. 

Now we focus on the anisotropic properties. In the 
main panel of Fig. 5, M L for the B || node (open cir- 
cles) and the B || antinode (filled circles) are shown. A 
small but clear difference in Ml between the two field 
orientations is seen. As we examine the data more care- 
fully, we find that the sign of the anisotropy changes 
between H c \ and H C 2. Namely, \Ml\ for the B | node 
is smaller near H C 2, while the tendency is reversed be- 
low field B* . This is more evident in Fig. 6, where the 
field-angle dependences of Ml are shown. Due to the in- 
herent fourfold anisotropy of d-wave pairing, Ml shows 
a fourfold oscillation as a function of the field angle 
Clearly, the sign of the M^-oscillation changes at field 
B*/jH orh ss 0.3 between H cl and H c2 . Thus, the sign 
change of Mi-oscillation can occur even in a simple d- 
wave superconductor. 
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Fig. 6. Field-angle dependences of Ml in a d-wave superconduc- 
tor. The parameters used are the same as those in Fig. 5. The 
data at different B are vertically shifted. 



In the following, we calculate the longitudinal magneti- 
zation Ml in a layered <i-wave superconductor under an 
in-plane field. We use the above formula for Ml and as- 
sume a triangular vortex lattice in the r'-frame. Since the 
calculation for an in-plane vortex state consumes much 
numerical resources, we use Simpson's rule for the spa- 
tial average as well as the Fermi surface average, and we 
set iV max = 6. 

Figure 5 shows the field dependence of M L in a <i-wave 
superconductor. Here, H orh = 1.037(|e|/7r)/27r^Q is the 
three-dimensional orbital limiting field in the isotropic 
case. We fix the temperature at T/T c — 0.5 and set 7 = 3 
for the uniaxial anisotropy parameter. From the figure, 
we can see the London behavior Ml oc h\{B C 2/B) in a 



3.3 Torque 

Next we discuss the transverse component of the mag- 
netization, or the torque. The transverse magnetization 
M T is given by M T = -B^ 1 ±(F /V), where x is the 
field angle shown in Fig. 3. In experiments, this quantity 
is mostly measured as a torque density, t/V = M x B 
or t/V — —j^(F/V). To obtain our numerical data of 

r, we first calculate the ^-dependence of F/V, then per- 
form a polynomial interpolation, and calculate r using 
the above relation. 

Figure 7 shows the field-angle dependences of torque 
r in a d-wave superconductor at T/T c = 0.5. As in the 
case of longitudinal magnetization, r shows a fourfold 
oscillation as a function of \- When the field is lowered 
from H c2 (~ 0.4757-f/orb in this case), the oscillation am- 
plitude first increases up to B/jH OI k ~ 0.3, as seen in 
Fig. 7(a), whereas it starts to decrease at lower fields, as 
seen in Fig. 7(b). 

To see the field dependence of the oscillation ampli- 
tude, we plot, in Fig. 8(a), the field dependence of r at 
a fixed angle x = tt/8- In the figure, we can see that the 
sign of the oscillation does not change between H c \ and 
H C 2- Although the accuracy of our numerical solution 
breaks down near B = or H = H c i, we believe that 
this behavior of r is not changed qualitatively by includ- 
ing the screening effect. Later, we will give a theoretical 
interpretation of the phenomena. It is worth noting that 
the behavior is in contrast to that in the case of longitu- 
dinal magnetization Ml, where the sign of the oscillation 
amplitude is reversed at field B* between H c i and H C 2. 
For comparison, we show, in Fig. 8(b), the corresponding 
field dependence of 5M L = M L (x = V 4 ) " m l(x = 0). 
It is interesting that near field B* where the sign of SMl 
is changed, the field dependence of the torque is also 
changed from a decreasing one to an increasing one. This 
coincidence is explained in the next subsection. 
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Fig. 7. Field-angle dependences of r in d-wave superconductor. 
The parameters used are the same as those in Fig. 5. Data in (a) 
high fields and (b) low fields are shown. 



o 
I 

> 



-0.001 



-0.002 



-0.003 



Q2D (y=3) 
d-wave 
T/Tc =0.5 




0.1 



0.2 0.3 

B/yHorb 



0.4 



0.5 



0.004 



y 
b 



to 



-0.004 




Q2D (y=3) 

d-wave 
T/Tc =0.5 



0.2 0.3 

B/yHorb 



0.5 
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used are the same as those in Fig. 5. 



3.4 Interpretation of results 

We discuss the interpretation of the results that the 
sign of Mi-oscillation changes, while the sign of r- 
oscillation does not. For this purpose, it_is useful to con- 
sider the anisotropy of the free energy F(x)- On the ba- 
sis of our numerical results, we schematically plot, in 
Fig. 9, the field dependences of F{x) in a d-wave super- 



-M L (B||node) 



F/V(B||node) 
F/V(B||antinode) 




Fig. 9. Schematic field dependences of F/V for B \\ antinode 
(solid line) and B || node (dashed line). 



conductor, corresponding to the two different field orien- 
tations (B || node and antinode). Here, the magnitude of 
anisotropy is enlarged for clarity. It is important to note 
that the two curves do not cross, and F/V is a convex 
function of B. This can be understood as follows. 

Firstly, it is important to note the relationship be- 
tween the anisotropies of coherence length and penetra- 
tion depth. 46, 50 Extending this knowledge, we obtain the 
anisotropy relation between iJ cl and H c2 : 21 

H cl ( X = 0) ff c2 (x = 7r/4) 



H cl ( X = n/4) " H c2 (x = 0) ' (45) 

This implies that if H c2 in a certain direction is larger 
than in other directions, then H c \ in that direction is 
smaller. This can be explained in an alternative way. Re- 
call that H cl = {\tik/V2k)H c and H c2 = (V2n)B c , 51 
where H c is the thermodynamic critical field. By defini- 
tion, H c does not depend on the field orientation, and we 
can ascribe the anisotropy to k. From these expressions, 
it is obvious that H cl and B c2 have opposite tendencies 
of anisotropies in the high-K case. 

Next we consider the anisotropy of F/V using the an- 
alytical results. 52 

-(H c2 -B) 2 /16ttk 2 (H^H c2 ) 
F(B = 0)/V + BH cl /4ir {H w H c i) 



F/V 



(46) 

Near H&, the anisotropy of F/V is determined by that 
of H c i, while near H c2 , the anisotropy of F/V is gov- 
erned by that of H c2 . Using eqs. (45) and (46), the free 
energy anisotropy shown in Fig. 9 can be confirmed an- 
alytically near H c \ and H c2 . In the intermediate field 
region, it is very difficult to imagine that the two curves 
F(x = 0) and F(x = 7r/4) cross at some field, i.e., the 
free energy becomes ^-independent. Otherwise, the d- 
wave superconductor appears effectively isotropic at this 
particular field. In this way, the free energy anisotropy 
shown in Fig. 9 can be understood. Moreover, we can 
show that F/V is a convex function of B if we correctly 
include the vortex interaction energy. 52 

Now let us first check the sign reversal of Mi- 
oscillation using Fig. 9. Ml is given by Ml — 
— -§g(F/V), i.e., the negative slope of the field depen- 
dence of F/V. Then it is obvious from Fig. 9 that Ml 
near H c2 is larger for the B \\ node, while near H c \, this 
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tendency is reversed. This explains the sign change of 
-oscillation. On the other hand, the in-plane torque 
t = —-j^F in the simple d-wave case is roughly esti- 
mated by comparing F(x) in the two characteristic di- 
rections, x = (B | antinodc) and \ = tt/A (B \\ 
node). Namely, we have the simple relation r <~ — SF = 
F(B || node) - F(B || antinode). As is seen in Fig. 9, 
the sign of this quantity is not changed between H c \ 
and H c2 . Thus we immediately conclude that the sign of 
r-oscillation is not changed. In the same way, we can 
see that — r ~ 5F is maximum near field B* where 
the difference in the slope of the two curves in Fig. 9, 
i.e., -5M L = %Q^{B || node) - %^{B || antinode), 
changes its sign. In this way, we can understand why near 
B* , where the sign of SMr, is changed, the field depen- 
dence of the torque is also changed from a decreasing one 
to an increasing one. 

We end this section with some comments on the effect 
of Fermi surface anisotropics. This effect can be studied 
by introducing parameter (3, that describes the in-plane 
anisotropy of the Fermi surfaces, in a similar manner as 
in ref. 21. We performed a calculation for when this Fermi 
surface anisotropy essentially competes with the gap 
anisotropy, and confirmed that the observed tendency 
of anisotropy is reversed at around |/3| = 1.0 (similar to 
the difference in Figs. 3(a) and 3(b) of ref. 21). As the 
reader may note, however, the above argument for the 
sign change of Ml -oscillation and the sign preservation 
of r-oscillation is independent of the source of anisotropy. 
The key is the use of eq. (45), i.e., the fact that H c \ and 
H C 2 have the opposite tendencies of anisotropy, and this 
relation holds irrespective of the source of anisotropy. 
Thus, even if the system is strongly affected by the Fermi 
surface anisotropy, the statement that the direction with 
maximal H c2 is always energetically stable in the whole 
field range from H c \ to H C 2 is valid. The importance 
of this fact is discussed in the next section. Repeating 
the same argument as in the previous paragraph, we can 
say that the sign change of M/,-oscillation and the sign 
preservation of r-oscillation are universal. 

4. Conclusion 

In this paper, we proposed a new method of obtain- 
ing the "full solution" of the quasiclassical Eilenberger 
equations for high-K superconductors, in addition to in- 
troducing a reasonable "approximate solution" to exam- 
ine thermodynamic quantities. The approximate solution 
discussed here is a natural extension of the one devel- 
oped by Pesch 26 and which has been used recently by 
several authors. 29-33 We compared the approximate so- 
lution with the full solution or with Pesch's solution, and 
showed that our approximate method gives an excellent 
result, at least with respect to the thermodynamic quan- 
tities. Recently, we have extended the approximate solu- 
tion to the case including Pauli paramagnetism. 53 

We then applied our approach to the description of the 
in-plane magnetic anisotropies in a layered d-wave super- 
conductor such as high-T c cuprate. We clarified that the 
sign of the longitudinal magnetization oscillation changes 
between H c \ and H c2 , while the sign of the torque oscilla- 



tion does not. Based on the analytical expressions of the 
free energy and the anisotropy relation between H c i and 
H c2 , we showed that there is a close relationship between 
Mi-oscillation and r-oscillation. 

Furthermore, we revealed that the basal-plane torque 
data always indicate the maximal i? C 2-direction. The im- 
plication of the present result is that the basal-plane B c2 
anisotropies can be studied through torque measurement 
at low fields, even in materials with strong fluctuations, 
such as high-T c cuprates and organic superconductors, 
where the true H C 2 cannot be determined experimen- 
tally. Moreover, we can detect the node position if we 
limit ourselves to a d-wave superconductor possessing 
a less anisotropic Fermi surface. For example, in light 
of our results, the torque data for YBCO 22 are consis- 
tent with d x 2 _ y 2 -pairing. In this respect, it will be in- 
teresting to measure the basal-plane torque in the d- 
wave superconductor k-(ET) 2 Cu(NCS)2, since the node 
position concluded from the earlier theoretical analy- 
s j g 54-56 an( j foe recent experiments 57 ' 58 are controver- 
sial. 59 Also of interest is the measurement of the basal- 
plane magnetic anisotropics in the d-wave superconduc- 
tor CeCoIn 5 , since there again is a controversy concern- 
ing the node position. 14 ' 60 
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